Genetic algorithms for the numerical solution of 
variational problems without analytic trial functions 
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A coding of functions that allows a genetic algorithm to minimize functionals without analytic 
trial functions is presented and implemented for solving numerically some instances of variational 
problems from physics. 



I. INTRODUCTION 

The genetic algorithm (GA)i2i£ has become popular as 
as a reliable computerized method for solving problems 
from a wide range of domains, such as function optimiza- 
tion, handling them even in nonlinear, multidimensional 
search spaces. A conventional GA is a stochastic search 
method inspired in natural selection and the darwinian 
idea of the survival of the fittest in which the possible 
solutions of a problem (metaphorically the points in the 
search space) are coded as fixed-length strings of char- 
acters of an alphabet (usually binary) that resemble the 
chromosomes of alive beings. A GA evolves a population 
of search "points" chosen randomly appliying iteratively 
on them operators of selection, crossover and mutation 
for creating new populations (generations). 

Selection consists in giving a proportionally bigger 
number of offspring to the fitter individuals so the char- 
acteristics that make them better prevail. The combina- 
tion of this characteristics for generating new individuals 
is achieved through crossover, that is the interchange of 
portions of the strings of characters of two individuals 
paired randomly, giving birth to two new individuals for 
the next generation. In its simplest form in a GA all 
individuals are removed (die) after reproduction. 

The last iterative step consists in making random 
changes to the strings of individuals chosen with a small 
probability, which is named mutation after the natural 
process that it resembles. After some generations the in- 
dividuals tend to concentrate around the fittest "points" 
in the search space, so it can be said that all of the pro- 
cess was a way of optimizing the function employed to 
determine the fitting. 

The predominant kind of optimization problems at- 
tacked with GAs have been those in which the strings 
code literally points in a multidimensional space, where 
each dimension represents a parameter or variable of in- 
terest. When the potential solutions of a problem are 
functions (as it is the case for variational problems^), 
not points, the most popular GA approach developed 
to date has been that of choosing a set of analytic trial 
functions and combining them in the fittest way. There 
are two main ways for doing so: weightcning them, case 
where a string of weights is an individual, or using ge- 
netic programming (GP)* where the trial functions and 
the mathematical operators needed for combining them 



are the alphabet that gives shape to each member of the 
population. 

In this paper a way to directly represent numerical 
functions as strings (individuals) of a GA is presented, 
followed by its successful implementation on some in- 
stances of variational problems from physics. 



II. ANGULAR GENES ARE NOT REAL GENES 

Lets take <3 as the alphabet chosen to code individuals. 
In the literature of GAs (without taking GP into account) 
two main alphabets are usually discussed, the so called 
binary genes © = 0, 1 and the real genes (5 = 5R. Even 
when there is not any special restriction on <3 in the 
definition of the GA, only the implicit warning that it 
must facilitate the heredity of the fittest characteristics 
for obtaining acceptable results, it is curious how the 
attention in the field has been biased toward the binary 
and real alphabets. One of the goals of this paper is to 
emphazise the importance of focusing attention in other 
alphabets, in the extra information that is possible to get 
from them, specifically in an angular one, which will be 
called from now on angular genes. Aren't angular genes 
just real genes? The distinction made is based in the 
commonly forgotten fact that angles are not numbers^ 
they are an entity by themselves. 

The angular genes code piecewise functions as a string 
a of the angles between each consecutive pair of linear 
segments. For any combination of angles it is possible 
to scale and rotate the collection to fit the initial and 
final desired values. Taking y\ — y{x\) and jjn — v(xn) 
as the initial and final values of the piecewise function 
Vh = y(xk) to be represented in the range (x\,Xn) with 
N — 1 linear segments and k = 2, 3, . . . , N, the coding is 
defined as follows: 

R = v 7 (vn - yi) 2 + (xn - xi) 2 



(3 = tan l [(y N -yi)/{x N -xi)] 
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III. EXAMPLES 
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The presented coding was used to solve instances of 
four well known variational problems from physics. In 
each case the population used had a size of 100 with N = 
101, probability of mutation of 0.05 and 500 generations. 
The number of runs made for all cases was ten. For 
the three first cases it was not needed the help of s, so 
Si = 1, but not for the last where < Sj < 1. The 
algorithm was written in MATLAB and implemented in 
a personal computer with Pentium(R)4 CPU, 2.4GHz, 
448 MB RAM. 



K I % — ± \ 

53 ( fl / r ) sin X>j-7 + /3 (2) 
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Where — a < otj < a and s is a string of real numbers 
< Si < 1 that codes the relative length of each linear 
segment and together with a forms an individual. Having 
J]J and 101 we can further define: 



Axi = ( R / r ) Si cos I atj - 7 + 

3 = 1 



A. Curve of shortest distance in the Euclidian 
plane 

The functional to minimize is 



J = 



.).: 



{dyf 



whose known solution^ is the straight line y = ax + b. 
For the case x± = y± = = Q, xpj = 1 with minimum 
J = 1 , the average solution found by the algorithm with 
a = 0.005-7T was J = 1 + 1.66 x 10~5 with standard 
deviation of 1.2 x 10~6. 



Ay* = ( R / r ) 
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V t = tan 
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y'l = tan(a i+ i)(l + y-y^ +1 )/Axi 



(3) 



(4) 



Taking 0^ = Q(xi) as the evaluation in Xi of the func- 
tion that minimizes the functional, for the best found 
individual we have: 



a i+i — tan 
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(5) 



Equation JSJ is a measure of the error of the best ap- 
proximation found that clarifies the influence of a proper 
choice of a according to the problem. If a is too small 
the error can be surely surpassed but if it is too big the 
search space grows. 

Ignoring the differences in a tenth of a and 1 the search 
space explored has a size of approximately lO^ 1 x 
20 JV - 1 . 



B. Curve of minimum revolution area 

Considering two paralel coaxial wire circles to be con- 
nected by a surface of minimum area that is generated by 
revolving a curve y(x) about the x-axis, the functional to 
minimize is 



J 



2Tty^/{dxf + {dyf 



whose known solution 4 is the catenoid y — cosh(ax + 
b)/a. For the case — x\ — xjq = 0.5, y\ = = 1 with 
minimum J = 5.9917, the average solution found by the 
algorithm with a = 0.0057T was J = 5.9919 with standard 
deviation of 1.4 x 10~5. 



C. Fermat's principle 

According to Fermat's principle light will follow the 
path y(x) for which 

xn,Vn 



J= I n{x,y) v / (dx) 2 + {dyf 

is minimum when n is the index of refraction. When 
n = e v the solution is y = ln(a/ cos(x + b)). For the 
case —xi = xm = 1)2/1= Vn — 1 with minimum J — 
4.5749, the average solution found by the algorithm with 
a = 0.017T was J = 4.5752 with standard deviation of 
5.2 x 10-5. 
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D. The energies of the hydrogen atom 

The hydrogen atomi is the quantum system made of 
a proton and an electron whose energies, without taking 
into account the degeneracies, can be found minimizing 
the functional 



E n — - 

c 



h 2 



n(n — l)h 2 
2^ 



dr 



with c = J °° u^dr, u n {r) — rR n , u n (0) = 0, i? 2 is the 
probability distribution for the radial location of the elec- 
tron, q is its charge, [i = m e m p / (m e + m p ) the reduced 
mass of the system, £o the permitivity of free space and h 
the Planck's constant divided by 2ir. The energies of the 
system are ruled by the equation E n = — 13.6052 eV/n 2 . 

The algorithm was used to find the three first ener- 
gies of the system. In this case the result found by each 
run depends strongly in a right choice of rjy such that 
u n (r N ) 2 w 0. For the ground state E 1 = -13.6052 eV 
the best found was about eight times the Bohr radius, 
the average solution found by the algorithm was E\ = 
-13.5987 eV with a standard deviation of 0.028 eV. For 
n = 2, E 2 = —3.4014 eV, the best choice for rjy made 
was about fifteen times the Bohr radius and the average 
solution found by the algorithm was E 2 = —3. 42467 eV 
with standard deviation of 0.005 eV. For n — 3, E3 = 
— 1.5117 eV, the best choice for made was about 
twenty five times the Bohr radius, reaching an average 



solution of E3 = — 1.5103 eV with standard deviation of 
0.001 eV. In the three cases a = 0.005tt. 



IV. CONCLUSIONS AND FUTURE WORK 

The examples shown were chosen with demostrative 
purposes. Better aproximations for specific cases can be 
reached increasing N and improving the choice of cr, with 
the extra computational effort it implies. Even thought 
that it was showed the efficiency of the coding to mini- 
mize the functionals presented it will be necessary the de- 
velopment of a theory of difficulty to give a more concise 
explanation of the kind of problems that could be hard to 
solve using it, like those already existent for binary genes 
like deceptioni£ and NK landscapes 9 ' . Another useful fu- 
ture development will be that of general ways of handling 
problems with constraints. An important potential ap- 
plication of the kind of genetic algorithm presented would 
be in those cases where there are not analytic solutions 
available, like in many quantum systems. 
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